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ABSTRACT 

We have observed the J = 3 — 2 transition of N2H + and N2D + to investigate 
the trend of deuterium fractionation with evolutionary stage in three selected 
regions in the Infrared Dark Cloud (IRDC) G28.34+0.06 with the Submillimeter 
Telescope (SMT) and the Submillimeter Array (SMA). A comprehensible en- 
hancement of roughly 3 orders of magnitude in deuterium fractionation over the 
local interstellar D/H ratio is observed in all sources. In particular, our sample 
of massive star-forming cores in G28.34+0.06 shows a moderate decreasing trend 
over a factor of 3 in the iV(N2D + )/iV(N2H + ) ratio with evolutionary stage, a 
behavior resembling what previously found in low-mass protostellar cores. This 
suggests a possible extension for the use of the iV(N2D + )/./V(N2H + ) ratio as an 
evolutionary tracer to high-mass protostellar candidates. In the most evolved 
core, MM1, the N 2 H + (3 — 2) emission appears to avoid the warm region traced 
by dust continuum emission and emission of 13 CO sublimated from grain man- 
tles, indicating an instant release of gas-phase CO. The majority of the N2H + and 
N 2 D + emission is associated with extended structures larger than 8" (~ 0.2 pc). 

Subject headings: ISM: abundances — ISM: clouds — stars: formation 

1. Introduction 

In the early evolutionary stages of star formation process, sequential depletion of molecu- 
lar species on grain mantles nurtures a peculiar low-temperature chemistry due to the removal 
of important gas-phase reactants, starting with sulfur-bearing species and followed by even 
volatile molecules such as CO (Bergin & Tafalla 2007). Besides H 2 , N 2 is thought to be least 
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affected in this condensation process and results in an enrichment of its daughter products, 
NH 3 and N 2 H + (Bergin et al. 2002). The removal of the gas-phase CO also promotes ion- 
molecular reactions and induces a sharp increase in the abundance of deuterated molecules 
in dense cores (Millar, Bennett, & Herbst 1989). Indeed, an enhancement of 2 — 3 orders of 
magnitude in the D/H ratio in star-forming cores (Crapsi et al. 2005; Fontani et al. 2006; 
Pillai et al. 2007) over the local interstellar value of 1.51 x 10 -5 (Oliveira et al. 2003) has been 
observed. In particular, the deuterium fractionation of N 2 H + , D iiac = iV(N 2 D + )/./V(N 2 H + ), 
in low-mass star-forming cores shows an increasing trend with dynamical age in the prestellar 
phase (Crapsi et al. 2005) but a decreasing trend in the protostellar phase (Emprechtinger 
et al. 2009). Chemical models anticipate D iiac to be affected by a few factors such as the 
kinetic temperature and, for ion species, the electron abundance (Roueff et al. 2005) as well 
as the gas-phase CO abundance (Aikawa et al. 2005). The abundance of gas-phase CO is 
expected to decline in the prestellar phase through molecular depletion onto grain surfaces 
but to rise up in the protostellar phase through sublimation of ice mantles as the envelope 
warms up. The correlation between deuterium fractionation and CO depletion factor has 
been recognized in a compiled sample of prestellar and protostellar cores (Crapsi et al. 2005; 
Emprechtinger et al. 2009). In a subsample of Taurus cores, a better correlation is found 
and leads to the speculation of external environment being influential to the evolution of a 
core (Crapsi et al. 2005). On the other hand, there is no evidence of a consistent behavior 
of -Df ra c in the case of high-mass protostellar candidates although a clear but less dramatic 
enhancement has been observed in a number of massive protostellar objects (Fontani et al. 
2006; Pillai et al. 2007). The environs of high-mass protostars may not sustain low temper- 
ature long enough to build up deuterated species as abundant as those of low-mass objects. 
In this study, we investigated a possible trend in the deuterium fractionation of N 2 H + with 
dynamical age using a sample of massive star-forming cores from one single IRDC to reduce 
the environmental fluctuations among the selected cores. 

Infrared dark clouds (IRDCs) were first discovered by Infrared Space Observatory (ISO) 
and the Midcourse Space Experiments (MSX) through silhouette against the bright, diffuse 
infrared background emission of the Galactic plane (Egan et al. 1998; Carey et al. 1998, 
2000; Hennebelle et al. 2001; Rathborne, Jackson, & Simon 2006). Because of their large 
mass (M > 1O 3 M ), low temperature (T < 20 K), and high density (n H2 > 10 5 cm~ 3 ), 
IRDCs have been proposed to be in the earliest stage of massive star formation. At a 
distance of 4.8 kpc, IRDC G28.34+0.06 (hereafter G28) is associated with roughly 10 3 M Q 
in the infrared absorption region and contains several dense cores in different evolutionary 
stages (Carey et al. 2000; Wang et al. 2008). We have selected three dense cores in G28 
to form an evolutionary sequence: starting with MM9 in an early stage of mass collection 
with one weak continuum source, followed by MM4 in a stage of mass fragmentation with at 
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least five continuum sources, and MM1 in a later stage with embedded massive protostellar 
objects with a total luminosity of 10 3 L Q . All the selected regions are associated with water 
masers, indicating star- forming activities (Wang et al. 2006). 

Our study is enabled by the Arizona Radio Observatory (ARO) Submillimeter Telescope 
(SMT) and the Submillimeter Array 4 (SMA). The SMT offers good sensitivity to detect 
weak line emission in extended structures while the SMA can preferentially image compact 
structures in the regions of interest. 

2. Observations and Data Reduction 

We observed the emissions of N 2 H+ (3 - 2) at 279.511780 GHz and N 2 D+ (3 - 2) at 
231.321864 GHz towards the three selected regions in G28 with both the SMT and SMA. 

2.1. SMT Observations 

Single-dish observations of N 2 H+ (3 - 2) and N 2 D+ (3 - 2) toward G28-MM1, MM4, 
and MM9 were carried out on 2008 April 13 with the 10-meter SMT on Mount Graham, 
Arizona. The observations were performed in the beam-switching mode. The pointing 
centers followed the peak positions given by Rathborne et al. (2006): (a, 5)(J2000) = 
(18 h 42 m 52.10 s , -3°59'45.0") for MM1, (a, <f)(J2000) = (18 h 42 m 50.70 s , -4°3'15.0") for MM4, 
and (a, 5)(J2000) = (18 h 42 m 46.70 s , -4°4'8.0") for MM9. The primary beam is about 26" 
for N 2 H + and 32" for N 2 D + . The spectral resolution is 1 MHz, corresponding to a velocity 
resolution of 1.07 kms" 1 for N 2 H + and 1.30 kins" 1 for N 2 D + . The temperature scale T£ 
was obtained using standard vane calibration, and the main beam temperature, T mb , was 
derived through T mb = T^/ry m b with a main beam efficiency ^ m b = 0.75. The respective rms 
noise level is about T mb = 20 and 10 mK for the N 2 H + and N 2 D + data, respectively. Data 
reduction was done with the CLASS package. 
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2.2. SMA Observations 

Observations with the SMA were carried out with seven antennas in the compact con- 
figuration on 2008 June 17 for N2H + and September 30 for N2D + . The system temperature 
varied from 150 to 350 K during the N 2 H + observation and from 100 to 180 K during the 
N 2 D + observation. The correlator was set to have a spectral resolution of 410 kHz, equiv- 
alent to 0.44 and 0.53 kms -1 for N2H + and N 2 D + , respectively. The projected baselines 
ranged from 11 — 118 kA and 11 — 59 kA, which are insensitive to structures larger than 
8" (Wilner & Welch 1994). The full- width at half power (FWHP) of the primary beam is 
roughly 44" for N 2 H+ and 53" for N 2 D+. 

The observing cycle comprised scans of 1751+096, MM1, MM4, MM9, and 1911-201, 
and was repeated every 27 minutes. The phase centers of the first observation run were 
the same as those of the SMT observations. Based on the results of the first observation 
run, we adjusted the phase center for the second observation run on the continuum peak 
in MM1 to (a, 5)(J2000) = (18 h 42 m 52.00 s , -3°59'53.00") and in MM9 to (a, 5)(J2000) 
= (18 h 42 m 46.46 s , — 4°4'15.1"). Data inspection, bandpass and flux calibrations, as well as 
temporal gain derivation were done within the IDL superset MIR. The flux scale was derived 
by observations of Uranus and is estimated to be accurate within 15%. 

Imaging was performed using MIRIAD package with natural weighting. To keep com- 
parison easy, we present all the maps with center positions consistent with the phase centers 
of the second observation run. For every source in each observation run, we used visibilities 
from both sidebands, each of a 2 GHz bandwidth, to generate a line-free continuum map. 
Channel maps were made with visibilities gridded into a velocity resolution of 1 kms -1 and 
resulted in rms noises, o = 0.31 and 0.20 K, equivalent to 110 and 82 mJybeam" 1 , with 
angular resolutions of 3.3" x 1.7" and 4.0" x 2.4" for N 2 H + and N 2 D + , respectively. 

3. Results and Discussion 

3.1. Deuterium Fractionation as an Evolutionary Probe 

The J = 3 — 2 transition of N 2 H + and N 2 D + were detected in all three sources with the 
SMT (Fig. [TJ. Since both transitions contain numerous hyperfine components, we fit each 
spectrum of every source with a model comprised of thirty-eight hyperfine components with 
updated line frequencies and spontaneous emission rates (Pagani, Daniel, & Dubernet 2009). 
For each individual source, all the hyperfine components of every J-level are assumed to be 
in thermal equilibrium at a single excitation temperature, T ex , adopted from the ammonia 



-5 - 



observations with an angular resolution of 40" (Pillai et al. 2006). The models are described 
by three more parameters: total column density, N, systemic velocity, "Ulsr,) full- width at half 
maximum (FWHM) as line width, Av. Model spectra are optimized with the minimization 
of the reduced \ 2 value, x 2 -, an d the results are listed in Table [TJ The N 2 H + spectrum 
towards MM1 appears to be doubly peaked, possibly affected by the presence of multiple 
sources (Zhang et al. 2009) and resulted in a significantly broader line width and a large 
X 2 - 

In cold clouds, a correction of the cosmic background temperature, T bg = 2.7 K, is 
necessary when extracting the optical depth information of an observed spectrum with 



t(v) = — In 



J(T cx )-J(T hg )\ 



(1) 



where T m b(f) is the main beam temperature of the spectra, and J(T) = (hv /k){e hv ' kT — 1). 
All the line emissions appear to have fairly small optical depths. Meanwhile, the op- 
timized model also provides an estimate of optical depth by integrating optical depths 
of all the hyperfine components. The maximum optical depth, r max , is less than 0.21 
(Table [T]). The emission of all the observed lines is optically thin. Over all, the fitted 
line widths, Av, of all the transitions are much broader than their thermal line width, 
A-Uth = a/8 In 2 kT cx / m N2 H+ — 0.16 kms -1 , suggesting a significant contribution from non- 
thermal motions. Our N2H + (3 — 2) spectra also have larger line widths when compared to 
the NH 3 (1, 1) spectra observed with the Very Large Array (VLA) and the Effelsberg 100m 
telescope (Pillai et al. 2007; Wang et al. 2008). Since the two transitions have similar upper 
state energies, E up = 26.8 K for N2H + and 23.4 K for NH3, but fairly different critical den- 
sities, n crit ~ 10 6 cm -3 for N 2 H + (Daniel et al. 2005; Pagani et al. 2009) and 10 3 cm -3 for 
NH 3 (Evans 1999), our N 2 H + observations tend to trace denser clumps possibly embedded 
in the inner region that is affected by star- forming activities as suggested by H2O masers 
(Wang et al. 2006). Observations of lower J transitions of N 2 H + such as the J = 1 — 
{E np = 4.5 K and n cr ; t ~ 10 5 cm -3 ) line are needed to verify this interpretation. 

The deuterium fractionation, -Df rac , is calculated for every source (Tabled]) and shows a 
significant enhancement of roughly 3 orders of magnitude higher than the average D/H ratio 
of 1.51 x 10~ 5 in the local interstellar medium (Oliveira et al. 2003). Such an enhancement 
has been observed in a large sample of high-mass star-forming cores (Fontani et al. 2006) 
as well as low-mass prestellar and protostellar cores (Crapsi et al. 2005; Roberts & Millar 
2007; Emprechtinger et al. 2009). In particular, we find a decreasing trend with evolutionary 
stage, from -Df rac = 0.051 in the younger MM9 core to -Df rac = 0.016 in the more evolved 
MM1 core, a change over a factor of roughly 3. Although a decreasing trend over a factor 
of roughly 8 has been identified in low-mass protostellar cores (Emprechtinger et al. 2009), 
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there was no evidence for such trend in a sample of massive star-forming cores associated 
with different molecular clouds (Fontani et al. 2006). Variations in thermal history, external 
environment, and initial chemical abundances across different molecular clouds may cause 
undesired confusion. Although our sample contains only three sources, the association in 
one single IRDC helps to minimize the environmental variations among the sources and to 
reveal the gentle drop in deuterium fractionation as the dusty envelope warms up. 



3.2. Emission Contained in Compact Structures 

Using the SM A, we further imaged the N 2 H + and N 2 D + emission to study more centrally 
concentrated structures. The N2H + emission shows clear detections in MM1 and MM4 as well 
as a marginal detection in MM9 (Fig. [2]). Since an interferometer is insensitive to structures 
larger than the scale corresponding to its shortest projected baseline, our SMA observations 
serve as a spatial filter to probe the fractional flux contained in clumps smaller than 8" 
(Wilner & Welch 1994). A comparison of the integrated intensity observed with the SMA, 
Wsma, to that with the SMT, Hsmt, can be made by convolving the SMA maps with the 
SMT beam after masking out regions below —3a, which are believed to be artifacts induced 
by the lack of short baselines. For every source, we compare the N 2 H + integrated intensity 
observed with the SMA to that with the SMT (Fig. [TJ and estimate the fraction of the 
SMA integrated intensity, Wsma/Wsmt 5 which indicates the state of a mass concentration 
process (Table [2]). The small values of Wsma/Wsut suggest that most of the N 2 H + emission 
is in structures larger than 8". An increasing trend from MM9 to MM1 is found and agrees 
with the presumed evolutionary stage of the central sources. On the other hand, N 2 D + is 
not detected in all regions at a level of 4a ~ 2.3 Kkms -1 , translating to an upper limit of 
iV(N 2 D+) = 2.0 x 10 12 cm" 2 at T ex = 13 K. Compared to the SMT integrated intensity 
in MM9, this detection limit sets a maximum fraction of 19% for emission coming from 
compact structures. Given the high critical density of n cr it ~ 10 6 cm -3 , the large fraction of 
the N 2 H + (3 — 2) emission missed by the SMA observations strongly suggests the presence of 
cold and dense gas in scales larger than 8" (~ 0.2 pc). Similar to N 2 H + , most of the N 2 D + 
emission is in extended structures. 

Continuum emissions at 284.2 and 236.5 GHz are detected in all sources with similar 
morphology. For the first time, a compact continuum source is detected in MM9 (Fig. H]). 
Since the visibility coverages of the two observation runs were quite different, we made 
continuum maps with visibilities of projected baselines within 13 — 57 kA, the range that all 
sources have in common at the two frequencies. If the dust emission is optically thin, the 
continuum flux density F u oc K u B u (Td), where k v = 0.006 (z//245 GHz)' 3 cm 2 g _1 (Shepherd 
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& Watson 2002) is the dust opacity at the observing frequency v, and B u (Td) is the Planck 
function at a dust temperature T d . The opacity spectral index, (3 = Alogfi^/Alogf, can 
be measured by comparing flux densities at two observing frequencies. In cold clouds, the 
condition hv ~ feT^ makes the Rayleigh- Jeans approximation inappropriate. We estimate 
the opacity spectral index in the selected regions with 

_ \og{FjF Ul ) + log[(e^/^ - l)/( e W^ i)] 
log(z/ 2 /z/i) 

and the core mass with M core = F U D 2 / K u B u (T d ), where D is the distance of the source. 
Assuming thermal equilibrium between gas and dust over similar spatial scales in the cal- 
culations of (3 and M core , we adopt the gas temperature derived from the VLA ammonia 
observations with a comparable angular resolution of 5" x 3" (Zhang et al. 2009) to be the 
dust temperature, T d . The results are listed in Table [2j Because of the limited frequency 
span of 48 GHz, we note that a systematic uncertainty of 15% in flux density measurements 
will produce a fairly large uncertainty of A/3 ~ 1.2 in our calculations. 

The dust opacity spectral index shows an increasing trend from MM9 to MM1 as the 
central sources evolve although the multiplicity in MM4 may complicate the interpretation 
of its averaged (3. Nevertheless, a smaller value of (3 is observed in MM9 with respect to 
MM1. At millimeter wavelengths, the value of (3 is a good probe for the size distribution of 
dust grains (Miyake & Nakagawa 1993), and a smaller (3 suggests that MM9 is surrounded 
by larger dust grains as a result of grain growth in high-density environment. On the other 
hand, the larger (3 in the more evolved region MM1 may be attributed to possible changes 
in the size distribution and chemical composition of the surrounding dust grains, which have 
been exposed to strong radiation fields generated by the associated massive YSOs. 



3.3. Deficiency of N 2 H + in a Warm Region 

In cold, dense environment, CO tends to freeze out onto dust grains while N2H + and 
NH 3 suffer less from depletion (Bergin et al. 2002; Tafalla et al. 2004). Since CO is the 
major destroyer of molecular ions, its removal from the gas phase results in a subsequent 
enrichment of N 2 H + (Aikawa et al. 2005). In later evolutionary stages when the internal 
heating of YSOs becomes important, the temperature will rise up and lead to sublimation 
of ice mantle, which returns volatile species such as CO back to the gas phase. The warmer 
temperature together with the reappearance of gas-phase CO can alter the competition 
among chemical reaction routes and destroy N 2 H + that has been produced during the cold 
early phase (Roueff et al. 2005; Aikawa et al. 2005). 
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In MM1, a deficiency of N2H + is observed in the location of the dust continuum emission 
as well as the 13 CO (2 — 1) emission (Fig. [3]) in a previous study (Zhang et al. 2009). 
Observationally, a centrally heated temperature structure from 16 to 30 K is derived in 
MM1 over spatial scales of roughly 1 to 0.1 pc (Pillai et al. 2006; Zhang et al. 2009). 
Theoretically, a significant drop in -Df rac is predicted across this temperature range (Roueff 
et al. 2005), and an instant release of CO is expected due to a significant drop of the CO 
sublimation timescale from 10 s yr at T<i ~ 12 K to 0.1 yr at ~ 20 K (Collings et al. 
2003). This hypothesized release of CO in the warm region can be traced in MM1 with the 
13 CO emission and dust continuum emission. A contrast in chemical composition can occur 
across the boundary between the cold, outer part and the warm, inner part that has been 
altered by the newly released gas-phase reactants. This chemical contrast has recently been 
observed in AFGL 5142 by comparing emissions of N 2 H + and NH 3 (Busquet et al. 2009). 



4. Summary 

Our main findings are summarized as follows: 

1. A moderate decreasing trend over a factor of 3 in deuterium fractionation of N 2 H + 
with evolutionary stage. Such trend extends the use of the iV(N2D + )/iV(N2H + ) ratio 
as an evolutionary tracer to high-mass protostellar candidates. 

2. An N 2 H + void resulting in part from a warmer temperature of the core and an instant 
release of CO from sublimated grain mantles in the warm region of MM1. 

3. A large fraction of the N 2 H + (3 — 2) flux missed by the interferometer suggesting the 
presence of cold and dense gas over a rather large area (> 0.2 pc) in G28. 

This research is supported by National Science Council of Taiwan through grants NSC 
97-2112-M-001-006-MY2 and NSC 97-2112-M-007-006-MY3. 
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Table 1. Results of SMT Spectral Fits 



Parameters 


MMl 


MM4 


MM9 


T cx (K) 




16.0 


16.6 


13.2 


N 2 H+ (3 


-2): 








N (10 12 


cm -2 ) 


10.12 ±0.05 


5.14 ±0.04 


1.86 ±0.05 


"^lsr (kms r ) 


77.59 ±0.02 


78.63 ±0.02 


79.82 ± 0.04 


Av (kms -1 ) 


6.14 ±0.04 


3.88 ±0.04 


3.42 ±0.11 


X 2 




122.9 


9.2 


1.5 


''max 




O.ZL 


U. io 


U.Uo 


N 2 D+ (3 


-2): 








N (10 11 


cm -2 ) 


1.7 ±0.2 


1.2 ±0.3 


1.0 ±0.2 


t; L SR (kms r ) 


78.0 ±0.2 


79.6 ±0.4 


81.1 ±0.2 


Av (kms -1 ) 


2.6 ±0.5 


3.7 ±0.9 


2.2 ±0.6 


X 2 




1.8 


0.9 


1.4 


'max 




0.007 


0.003 


0.005 






0.017 ±0.002 


0.024 ±0.005 


0.052 ±0.011 




Table 2. Comparison of Three Sources 




Parameters 


MMl MM4 


MM9 




N 2 H+ (3-2): 








WsMA 


(Kkms -1 ) 


1.12 0.33 


0.05 




W^SMT 


(Kkms- 1 ) 


12.93 7.17 


2.23 




W SMA /W SMT (%) 


8.6 4.6 


2.4 




Continuum emission: 








T 

- 1 ex 




30.0 13.0 


13.7 




^236.5 GHz (Jy) 


0.59 0.15 


0.03 




^284.2 GHz (Jy) 


1.24 0.28 


0.04 




P 




2.3 2.0 


0.9 




■^core 




277 209 


36 
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Fig. 1.— (a) through (c) SMT N 2 H+ (3-2) spectra (histogram) towards MM1, MM4, 
and MM9, respectively, superposed with the thirty-eight hyperfine component model-fit 
spectra (solid curves) assuming a single excitation temperature. The frequency of each 
individual hyperfine component is labelled with a short bar on the top. The SMA spectra 
after convolving with the SMT beam are plotted with a displacement of —0.1 K (thick 
histogram), (c) through (f) SMT N 2 D + (3 — 2) spectra (histogram) towards MM1, MM4, 
and MM9, respectively, superposed with the thirty-eight hyperfine component model- fit 
spectra (solid curves). 
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Fig. 2. — SMA N2H + (3 — 2) integrated intensity [contours) overlaid on the 284.2 GHz line- 
free continuum map towards MM1, MM4, and MM9. Black crosses label the continuum peaks 
reported by Zhang et al. (2009) in MM1 and MM4 and by this work in MM9. Contour levels 
correspond to (-10, -8, . . . , -2, 2, 4, . . . , 22) a, where a = 1.10,0.93, and 0.88 Kkms' 1 for 
MM1, MM4, and MM9, respectively. The synthesized beam shown is 3.3" x 1.7" (P. A. = 49°) 
for the N2H + maps. The continuum map has an intensity range of —48 to 738 mJybeam -1 
for MM1, —14 to 56 mJybeam -1 for MM4, and —13 to 46 mJybeam -1 for MM9 with beam 
sizes of 4.0" x 2.5" (P. A. = 31°). 



- 14 - 



7 



G28 



-MMf - 



o 
o 
o 

CM 

O 
<D 
O 



50 



-04 00 00 



1 8"42 m 53r0 



52?5 52f0 

RA (J2000) 



5175 



5170 



Fig. 3. — SMA N2H + (3 — 2) integrated intensity (red contours) and 13 C0 (2 — 1) integrated 
intensity (white contours) overlaid with the line- free continuum map (colors) of MM1. Con- 
tour levels correspond to (—10, —6, —2, 2, 6, ... , 22) a, where a = 1.10 and 1.12 Kkms' 1 for 
N2H + and CO, respectively. Black crosses indicate the continuum peak positions (Zhang 
et al. 2009). The synthesized beam size is 4.0" x 2.8" (P.A. = 44°) for CO. 



